Задание 2

plot_ly(data = data_raw_0[(data_raw_0$`Sucide Rate` != 0) & (data_raw_0$`Urban population` != 0),],
        x     = ~ `Sucide Rate`,
        y     = ~ `Urban population`,
        color = ~ continent,
        marker = list(
        size = 10
      )
)   %>%
  layout(
    title = 'Отношение Sucide Rate к Urban population',
    yaxis = list(title = 'Urban population',
                 zeroline = FALSE),  # Уберём выделения нулевых осей по y
    xaxis = list(title = 'Sucide Rate',
                 zeroline = FALSE)) # Уберём выделения нулевых осей по y

Задание 3

data_raw_0 %>% 
  filter(continent == "Africa" |  continent == "Americas") %>%
  group_by(continent) %>%
  get_summary_stats(`Life expectancy`, type = "mean_sd")%>%
  ungroup()
## # A tibble: 2 × 5
##   continent variable            n  mean    sd
##   <fct>     <fct>           <dbl> <dbl> <dbl>
## 1 Africa    Life expectancy    52  66.1  5.97
## 2 Americas  Life expectancy    38  78.6  3.68
data_raw_filtred <- data_raw_0 %>% filter(continent == "Africa" |  continent == "Americas") 
ggqqplot(data_raw_filtred[data_raw_filtred$`Life expectancy`,], 
         x = "Life expectancy", facet.by = "continent")

Так как данные распределены “не нормально” предлагаю использовать параметрический тест, например U-тест Манна-Уитни

library(rstatix)
# Создадим два датафрейма с значениями Life expectancy для разных континентов

group_af <- data_raw_0 %>% filter(continent == "Africa")   %>% select(`Life expectancy`) %>% pull(`Life expectancy`)
group_am <- data_raw_0 %>% filter(continent == "Americas") %>% select(`Life expectancy`) %>% pull(`Life expectancy`)

# Проведем тест Манна-Уитни

p_valie_Utest <- wilcox.test(group_af, group_am)
p_valie_Utest
## 
##  Wilcoxon rank sum test with continuity correction
## 
## data:  group_af and group_am
## W = 107, p-value = 6.342e-13
## alternative hypothesis: true location shift is not equal to 0
wilcox_test_result <- wilcox.test(`Life expectancy` ~ continent, data = data_raw_filtred) 
 
# Визуализация результатов 
ggplot(data_raw_filtred, aes(x = continent, y = `Life expectancy`, fill = continent)) + 
  geom_boxplot() + 
  geom_dotplot(binaxis = "y", stackdir = "center", position = "dodge") + 
  stat_compare_means(test = "wilcox.test", label = "p.format", comparisons = list(c("Africa", "Americas"))) + 
  labs(title = "Comparison of Life Expectancy between Africa and Americas") + 
  theme_minimal()

Задание 4

# Выбираем все числовые колонки, кроме 'Year' 
data_raw_0 <- data_raw_0 %>%
  mutate(continent_num = case_when(continent == "Africa" ~ 1,
                                   continent == "Americas" ~ 2,
                                   TRUE ~ 0))

numerical_data <- data_raw_0 %>% select(where(is.numeric)) %>% filter(continent_num > 0)
numerical_data <- subset(numerical_data, select = -c(Year))
 
# Проводим корреляционный анализ 
correlation_matrix <- cor(numerical_data) 
 
# Визуализация корреляций - используем heatmap 

 
melted_correlation <- melt(correlation_matrix) 
 
ggplot(melted_correlation, aes(x = Var1, y = Var2, fill = value)) + 
  geom_tile(color = "white") + 
  scale_fill_gradient2(low = "blue", high = "red", mid = "white", midpoint = 0, limit = c(-1, 1), space = "Lab", name="Correlation") + 
  theme_minimal() + 
  theme(axis.text.x = element_text(angle = 45, vjust = 1, size = 8, hjust = 1)) + 
  labs(title = "Correlation Heatmap") 

corrplot(correlation_matrix, method = 'ellipse', order = 'AOE', tl.cex = 0.5, tl.srt = 45)

Задание 5

numerical_data_scaled <- scale(numerical_data)
dist <- dist(numerical_data_scaled, 
                        method = "euclidean"
                        )
as.matrix(dist)[1:6,1:6]
##          1        2        3        4        5        6
## 1 0.000000 6.074379 5.025993 5.170680 5.884635 3.816022
## 2 6.074379 0.000000 8.099954 7.772997 6.845781 6.879864
## 3 5.025993 8.099954 0.000000 5.590110 4.748777 4.197357
## 4 5.170680 7.772997 5.590110 0.000000 4.719991 3.555197
## 5 5.884635 6.845781 4.748777 4.719991 0.000000 4.741891
## 6 3.816022 6.879864 4.197357 3.555197 4.741891 0.000000
clear_hc <- hclust(d = dist, method = "ward.D2")

fviz_dend(clear_hc, 
          cex = 0.1)
## Warning: The `<scale>` argument of `guides()` cannot be `FALSE`. Use "none" instead as
## of ggplot2 3.3.4.
## ℹ The deprecated feature was likely used in the factoextra package.
##   Please report the issue at <https://github.com/kassambara/factoextra/issues>.
## This warning is displayed once every 8 hours.
## Call `lifecycle::last_lifecycle_warnings()` to see where this warning was
## generated.

Задание 6

# Создаем heatmap с дендрограммой 
pheatmap(numerical_data_scaled, 
         show_rownames = FALSE, 
         clustering_distance_rows = dist,
         clustering_method = "ward.D2", 
         cutree_rows = 3,
         cutree_cols = length(colnames(numerical_data_scaled)),
         angle_col = 45, 
         main = "Dendrograms for clustering rows and columns with heatmap")

Задание 7

Делаем PCA:

full.pca <- prcomp(numerical_data_scaled, 
                        scale = F)

summary(full.pca)
## Importance of components:
##                           PC1    PC2     PC3     PC4     PC5     PC6     PC7
## Standard deviation     2.7815 1.6834 1.40851 1.36879 1.10759 0.98271 0.84961
## Proportion of Variance 0.3868 0.1417 0.09919 0.09368 0.06134 0.04829 0.03609
## Cumulative Proportion  0.3868 0.5285 0.62774 0.72142 0.78275 0.83104 0.86713
##                            PC8     PC9    PC10    PC11    PC12    PC13    PC14
## Standard deviation     0.74926 0.70909 0.65772 0.55364 0.47697 0.43659 0.41107
## Proportion of Variance 0.02807 0.02514 0.02163 0.01533 0.01138 0.00953 0.00845
## Cumulative Proportion  0.89520 0.92034 0.94197 0.95730 0.96867 0.97820 0.98665
##                          PC15    PC16    PC17    PC18    PC19      PC20
## Standard deviation     0.3317 0.27505 0.22546 0.16374 0.06028 2.257e-16
## Proportion of Variance 0.0055 0.00378 0.00254 0.00134 0.00018 0.000e+00
## Cumulative Proportion  0.9921 0.99594 0.99848 0.99982 1.00000 1.000e+00
fviz_eig(full.pca, addlabels = T, ylim = c(0, 50), ncp = 19)

fviz_pca_var(full.pca, 
             select.var = list(contrib = 5), # Задаём число здесь 
             col.var = "contrib")

Наибольший вклад в “коррелированность” данных вносят Life expectancy, Infant Mortality, HepB3 Immunization, DPT Immunization, Measles Immunization

fviz_contrib(full.pca, choice = "var", axes = 1, top = 30)

fviz_contrib(full.pca, choice = "var", axes = 2, top = 24)

fviz_contrib(full.pca, choice = "var", axes = 3, top = 24) 

fviz_contrib(full.pca, choice = "var", axes = 4, top = 24) 

fviz_contrib(full.pca, choice = "var", axes = 5, top = 24)

Рассмотрим Dim-1. В данных нет ярковыраженных переменных, которые бы брали на себя наибольших “вклад” для первой главной компоненты однако можно выделить группу переменных, которые отвечают за приблизительно 40% такого вклада:

  1. Life expectancy

  2. Basic sanitation services

  3. Infant Mortality

  4. Clean fuels and cooking technologies

Задание 8, 9

# Визуализируем с группировкой по continent (для этого переменную нужно сделать фактором)
ggbiplot(full.pca, 
         scale=0, 
         groups = as.factor(numerical_data$continent_num), 
         ellipse = T,
         alpha = 0.2) +
  theme_minimal()

Значащие группы 1 и 2 – это Африка и Америки соответственно.

Если мы доверяем построиным линиям то можно сформулировать следубщие гипотрезы:

  1. Страны с большим Sucide Rate более вероятно будут относиться к континенту африка

  2. Страны с больней Life expectancy более вероятно относятся к континантам Америки

  3. Rural population больше в странах Африки, Urban population больше в странах относящихся к Америке

  4. Показатель Clean fuels and cooking technologies и Basic sanitation services более выражены в странах Америки

и другие

Замечание не получилось перевести в plotly :(
#install.packages(c("irlba", "Matrix"))

library(tidymodels)
library(embed)

umap_prep <- recipe(~., data = numerical_data) %>% # "техническая" строка, нужная для работы фреймворка tidymodels
  step_normalize(all_predictors()) %>% # нормируем все колонки
  step_umap(all_predictors()) %>%  # проводим в UMAP. Используем стандартные настройки. Чтобы менять ключевой параметр (neighbors), нужно больше погружаться в машинное обучение
  prep() %>%  # "техническая" строка, нужная для работы фреймворка tidymodels. Мы выполняем все степы выше 
  juice() # Финальная строка - приводим результаты UMAP к стандартизированному датасету

Задание 10

umap_prep %>%
  ggplot(aes(UMAP1, UMAP2)) + #  # можно добавить раскраску 
  geom_point(aes(color = as.factor(numerical_data$continent_num),
                 shape = numerical_data$diabetes_ch), 
             alpha = 0.7, size = 2) +
  labs(color = NULL)+
  ggtitle("Отображение данных в зависимости от Континента")

Сложно интерпритировать UMAP график. Однако если сравнивать его с PCA графиком, то можно заметить что:

PCA: Основан на линейных преобразованиях и стремится сохранить максимальную дисперсию данных в новых компонентах (главных компонентах).

UMAP: Основан на нелинейных преобразованиях, стремится сохранить локальные структуры данных, предоставляя более низкоразмерное представление, сохраняя при этом глобальные структуры.

PCA: Сохраняет глобальные расстояния, что может привести к сжатию или искажению локальных структур данных.

UMAP: Стремится сохранять локальные расстояния, что делает его лучшим для визуализации кластеров и плотных областей в данных.

Задание 11

numerical_data_scaled_drop_1_5   <- numerical_data_scaled[, -c(1:5)]
numerical_data_scaled_drop_6_10  <- numerical_data_scaled[, -c(6:10)]
numerical_data_scaled_drop_11_15 <- numerical_data_scaled[, -c(11:15)]

full_drop_5.pca <- prcomp(numerical_data_scaled_drop_1_5, 
                        scale = F)

summary(full_drop_5.pca)
## Importance of components:
##                           PC1    PC2    PC3     PC4     PC5     PC6     PC7
## Standard deviation     2.4452 1.6408 1.2933 1.04801 0.86905 0.84416 0.77616
## Proportion of Variance 0.3986 0.1795 0.1115 0.07322 0.05035 0.04751 0.04016
## Cumulative Proportion  0.3986 0.5781 0.6896 0.76281 0.81316 0.86067 0.90083
##                            PC8     PC9    PC10    PC11    PC12    PC13    PC14
## Standard deviation     0.68995 0.58551 0.50801 0.45351 0.31760 0.27616 0.16689
## Proportion of Variance 0.03174 0.02285 0.01721 0.01371 0.00672 0.00508 0.00186
## Cumulative Proportion  0.93256 0.95542 0.97262 0.98633 0.99306 0.99814 1.00000
##                             PC15
## Standard deviation     2.294e-16
## Proportion of Variance 0.000e+00
## Cumulative Proportion  1.000e+00
full_drop_10.pca <- prcomp(numerical_data_scaled_drop_6_10, 
                        scale = F)

summary(full_drop_10.pca)
## Importance of components:
##                           PC1    PC2    PC3    PC4     PC5     PC6     PC7
## Standard deviation     2.3863 1.4573 1.3644 1.2733 1.01672 0.91987 0.80672
## Proportion of Variance 0.3796 0.1416 0.1241 0.1081 0.06891 0.05641 0.04339
## Cumulative Proportion  0.3796 0.5212 0.6453 0.7534 0.82231 0.87872 0.92211
##                            PC8     PC9    PC10    PC11    PC12    PC13    PC14
## Standard deviation     0.65115 0.48608 0.43437 0.40561 0.29331 0.25535 0.06080
## Proportion of Variance 0.02827 0.01575 0.01258 0.01097 0.00574 0.00435 0.00025
## Cumulative Proportion  0.95037 0.96612 0.97870 0.98967 0.99541 0.99975 1.00000
##                             PC15
## Standard deviation     2.879e-16
## Proportion of Variance 0.000e+00
## Cumulative Proportion  1.000e+00
full_drop_15.pca <- prcomp(numerical_data_scaled_drop_11_15, 
                        scale = F)

summary(full_drop_15.pca)
## Importance of components:
##                          PC1   PC2    PC3     PC4     PC5     PC6     PC7
## Standard deviation     2.456 1.464 1.3556 1.15139 1.02849 0.82952 0.72496
## Proportion of Variance 0.402 0.143 0.1225 0.08838 0.07052 0.04587 0.03504
## Cumulative Proportion  0.402 0.545 0.6675 0.75584 0.82636 0.87223 0.90727
##                            PC8     PC9    PC10    PC11    PC12    PC13    PC14
## Standard deviation     0.65098 0.60165 0.45755 0.40791 0.39998 0.25634 0.06107
## Proportion of Variance 0.02825 0.02413 0.01396 0.01109 0.01067 0.00438 0.00025
## Cumulative Proportion  0.93552 0.95966 0.97361 0.98471 0.99537 0.99975 1.00000
##                             PC15
## Standard deviation     2.384e-16
## Proportion of Variance 0.000e+00
## Cumulative Proportion  1.000e+00

Действительно Cumulative Proportion отличается в зависимости от выбранных переменных. Отличия не большие +- 5 процентных пунктов в моем случае.

Задание 12

Давайте самостоятельно увидим, что снижение размерности – это группа методов, славящаяся своей неустойчивостью. Создайте две дамми-колонки о том: (1) принадлежит ли страна к африканскому континенту, (2) Океании. Проведите PCA вместе с ними, постройте биплоты. Проинтерпрейтируйте результат. Объясните, почему добавление дамми-колонок не совсем корректно в случае PCA нашего типа.

data_raw_dummy <- data_raw_0 %>%
  mutate(is_africa  = ifelse(continent== "Africa",   1, 0),
         is_americas = ifelse(continent == "Americas", 1, 0)) %>%
  filter(continent_num > 0)
  
  
  data_raw_dummy <- subset(data_raw_dummy, select = -c(Year)) %>%
  select(where(is.numeric))
data_raw_dummy_scaled <- scale(data_raw_dummy)
dist <- dist(data_raw_dummy_scaled, 
                        method = "euclidean"
                        )
as.matrix(dist)[1:6,1:6]
##          1        2        3        4        5        6
## 1 0.000000 6.074379 5.776495 5.902815 6.537294 4.761230
## 2 6.074379 0.000000 8.585834 8.278090 7.414311 7.445792
## 3 5.776495 8.585834 0.000000 5.590110 4.748777 4.197357
## 4 5.902815 8.278090 5.590110 0.000000 4.719991 3.555197
## 5 6.537294 7.414311 4.748777 4.719991 0.000000 4.741891
## 6 4.761230 7.445792 4.197357 3.555197 4.741891 0.000000
dummy.pca <- prcomp(data_raw_dummy_scaled, 
                        scale = F)
summary(dummy.pca)
## Importance of components:
##                           PC1    PC2     PC3     PC4     PC5     PC6     PC7
## Standard deviation     3.0154 1.6938 1.41620 1.38738 1.22256 0.98898 0.88165
## Proportion of Variance 0.4133 0.1304 0.09116 0.08749 0.06794 0.04446 0.03533
## Cumulative Proportion  0.4133 0.5437 0.63488 0.72237 0.79031 0.83477 0.87010
##                            PC8     PC9    PC10    PC11    PC12    PC13    PC14
## Standard deviation     0.79926 0.71470 0.65792 0.57864 0.55150 0.43942 0.41479
## Proportion of Variance 0.02904 0.02322 0.01968 0.01522 0.01383 0.00878 0.00782
## Cumulative Proportion  0.89914 0.92236 0.94203 0.95725 0.97108 0.97985 0.98767
##                          PC15    PC16    PC17    PC18    PC19      PC20
## Standard deviation     0.3349 0.27870 0.22559 0.16379 0.06028 2.534e-16
## Proportion of Variance 0.0051 0.00353 0.00231 0.00122 0.00017 0.000e+00
## Cumulative Proportion  0.9928 0.99630 0.99862 0.99983 1.00000 1.000e+00
##                            PC21      PC22
## Standard deviation     1.47e-16 3.329e-32
## Proportion of Variance 0.00e+00 0.000e+00
## Cumulative Proportion  1.00e+00 1.000e+00
fviz_eig(dummy.pca, addlabels = T, ylim = c(0, 50), ncp = 19)

первая и последующие переменные стали “объяснять” меньший процент корреляции данных

fviz_pca_var(dummy.pca, 
             select.var = list(contrib = 5), # Задаём число здесь 
             col.var = "contrib")

первые 5 переменных не изменились но направление изменилось

# Визуализируем с группировкой по continent (для этого переменную нужно сделать фактором)
ggbiplot(dummy.pca, 
         scale=0, 
         groups = as.factor(data_raw_dummy$continent_num), 
         ellipse = T,
         alpha = 0.2) +
  theme_minimal()

Добавив бинарные дамми колонки которые характеризуют наблюдения в зависимости от пренадлежности к материку: Африка или Америка(Я не стал брать Океанию для консистенси с прошлыми графиками) похоже что таким образом мы “отделили” наблюдения одного материка от другого, но таким образом потеряли визуальное подтверждение “схожести” некоторых наблюдений (пересечение овалов на предыдущем ggbiplot)